Robustness of a Cellular Automata Model for 

the HIV Infection l 



oo 

P. H. Figueiredo, a '* S. Coutinho, b R. M. Zorzenon dos Santos 



in 



a Departamento de Fisica, Universidade Federal Rural de Pernambuco, Dois 
Irmdos, CEP 52171-900, Recife, Pernambuco, Brazil. 

h Laboratorio de Fisica Teorica e Computacional, Departamento de Fisica, 
Universidade Federal de Pernambuco, Cidade Universitdria, CEP 50670-901, 
i-5 '. Recife, Pernambuco, Brazil. 






C*2 



£ 



Abstract 

An investigation was conducted to study the robustness of the results obtained 
from the cellular automata model which describes the spread of the HIV infection 
within lymphoid tissues [1]. The analysis focussed on the dynamic behavior of the 
model when defined in lattices with different symmetries and dimensionalities. The 
results illustrated that the three-phase dynamics of the planar models suffered mi- 
pq ■ nor changes in relation to lattice symmetry variations and, while differences were 

1^ \ observed regarding dimensionality changes, qualitative behavior was preserved. A 

C further investigation was conducted into primary infection and sensitiveness of the 

latency period to variations of the model's stochastic parameters over wide ranging 
values. The variables characterizing primary infection and the latency period exhib- 
ited power-law behavior when the stochastic parameters varied over a few orders 
^ ■ of magnitude. The power-law exponents were approximately the same when lattice 

symmetry varied, but there was a significant variation when dimensionality changed 
from two to three. The dynamics of the three-dimensional model was also shown to 
be insensitive to variations of the deterministic parameters related to cell resistance 
to the infection, and the necessary time lag to mount the specific immune response 
to HIV variants. The robustness of the model demonstrated in this work reinforce 
that its basic hypothesis are consistent with the three-stage dynamic of the HIV 
infection observed in patients. 
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1 Introduction 



Since the first cases of Acquired Immunodeficiency Syndrome (AIDS) were 
notified during the early 1980s, considerable advances have been made in un- 
derstanding the dynamics of the infection generated by the Human Deficiency 
Virus (HIV) . What differentiates HIV from many other viruses is its ability to 
remain in the system after primary infection. Although the infectious process 
varies from patient to patient, a typical infection pattern has been observed 
amongst patients [213] and exhibits three phases: primary infection, the la- 
tency period and the onset of AIDS. Primary infection is characterized by a 
broad viral dissemination, which declines markedly (without eliminating the 
virus)in weeks to months after the emergence of the first HIV-specific immune 
response jl] . This rapid decrease in the plasma viremia titer is followed by the 
clinical latency period, which may vary from 2 to 10 years or more, with a very 
low viral burden. During the clinical latency period, although many patients 
are usually asymptomatic, all suffer a gradual performance deterioration of 
the immune systems, manifested by the decrease of CDA + T cell counts j2] 
to very low concentrations. The onset of AIDS is associated to achieving very 
low thresholds of CDA + T cell counts corresponding to 20-35% of the regular 
values of a normal individual. In the absence of treatment, the outcome of this 
stage is death by opportunistic diseases. 

The main target of the HIV is the CD4 + T cell. The viral envelop protein 
binds to some specific cell receptor, thus enabling the virus to deposit its ge- 
netic material into the cell. Once inside, it uses the host cell machinery to make 
copies of its viral DNA in the same manner as other retroviruses [5]. When 
the host cells are able to produce both the DNA and the viral protein envelop, 
they are able to release new viral particles into the extra-cellular environment. 
Replication of HIV is extremely fast and mutation rates are very high. It is 
estimated that during the HIV reproduction process, two generations are pro- 
duced per week with one mutation per generation [3] . Part of the new strains 
would be recognized by previously developed immune responses, but part of 
the variants would not be recognized and would activate new specific immune 
responses to them. Since the virus is not eliminated, the process of contin- 
uously activating a new specific immune response is maintained indefinitely, 
thus placing stress on the immune system. The ultimate manifestation of the 
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effects of this stress is a decline in T cell counts. Several mechanisms have 
been proposed to explain the decline of T-cell counts [2|6] the aggregation of 
infected and non-infected cells in the lymphoid tissues (syncytia formation), 
the destruction of memory-T cells, programmed cell death associated to lym- 
phocyte activation instead of cell proliferation, and high viral replication and 
mutation rates. However, despite all efforts to determine the cause of this 
decline during recent decades, the question still remains open-ended. 

Over the past two decades many mathematical models based on differential 
equation approaches have been proposed in order to study the different as- 
pects involved in the development of the HIV infection [7j. Most of these 
models however, describe the evolution of the virus and cell populations in a 
compartmentalized manner, and do not take into account the spatial localiza- 
tion necessary to mount the specific immune responses [5]. Although this kind 
of approach has contributed to the understanding of different aspects of HIV 
infection dynamics, it has not been able to reproduce the entire three-stage 
dynamics observed in infected patients. Depending on the adopted approach, 
the model only adequately describes either the early (stochastic models) or 
later stages (deterministic models) of the dynamics [7] . There is also a chance 
it may describe the two time scales using different sets of parameters for dif- 
ferent stages [8] . Discrete approaches based on binary cellular automata have 
also been proposed [9] to describe the interactions amongst CDA + and CD8 + 
T cells, macrophages and viral particles. These models reproduce the different 
attractors of the dynamics [TU] or the early stages of the infection [§] but not 
its entire course. Nowadays, spatially structured models are being recognized 
as a step forward in understanding the dynamics of virus infection processes, 
in which the influence of local dispersal of the virus and virus-target cells 
is relevant for disease persistence in vivo. This is the case of the HIV infec- 
tion [TTlfT2j . 

Recently, two of the authors have proposed a two-dimensional stochastic cellu- 
lar automata (CA) model in order to describe the spread of the HIV infection 
PQ, and which takes into account the spatial localization of target cells (T 
cells) that occurs in the lymphoid tissues to mount the specific immune re- 
sponse. In the case of the HIV infection it is this localization that contributes 
to the spread of infection throughout the cells. The model illustrates that 
the combination of the timescale involved in the immune response of any 
(healthy) individual with fast replication and high mutation rates of the HIV, 
together with the spatial localization generated in the lymphoid tissue, may 
go some way to explaining the three-stage dynamics observed in experimental 
findings [2|5] . By reproducing the two timescales observed in the experimen- 
tal data (concerning T cell counts and viremia titer) the model permits the 
possible mechanisms underlying the dynamics of this infectious process to be 
investigated. The results suggest that the slow timescale is associated to the 
agglutination of infected cells in structures, which may compromise the entire 



tissue and trap healthy cells, leading to the onset of AIDS. All these struc- 
tures may be associated to syncytia formation. Syncytia are structures, which 
are formed when infected blood is mixed with healthy blood, in in vitro ex- 
periments. Therefore, results suggest that what is observed in vitro may also 
occur in vivo. 

The cellular automata model proposed in reference [1] is defined on a two- 
dimensional square lattice so as to mimic the structure of lymphoid tissues 
where cell interactions take place. In the case of HIV infection, lymph nodes 
play a major role within lymphoid tissues. The lymph nodes are small organs 
exhibiting a bean shape. From electromicrography images [13] the region where 
the interactions between the different cells and the viruses take place, has 
a fractal structure similar to that of a sponge. The cells take in the order 
of hours to cross the lymphoid tissue and find the right place to interact 
within the porous structure. In the model, the range of interaction between 
the cells is approximated by a surface using a square lattice. Despite the 
approximation and remarkable accuracy of the model in describing the two 
timescales and the three-stage dynamics, the question concerning the relevance 
of the lattice geometry as well as its dimensionality (2-D, 3-D or fractal) 
remains an important issue to be investigated. This issue was taken up and 
addressed by Figueiredo [H] whose results are for the first time reported in 
this paper. Recently, Omerond [15] conducted a similar investigation into how 
different tilings of square and cubic lattices, and definitions of the nearest 
neighbors qualitatively affect the dynamics of the system, without exploiting 
the robustness. However, while being a more qualitative investigation, it was 
less detailed than the current study. 

Below, further results will also be presented concerning the robustness of the 
dynamics of the model when varying the stochastic parameters of the model. 
To be more precise, studies were undertaken into the characteristic changes 
of the primary infection, due to variations of the initial HIV concentration 
for both the two and three-dimensional models. Moreover, an analysis was 
carried out of the time-scaling behavior of the clinical latency period based 
on the probability that newly infected cells enter the system throughout time, 
by different lattices. In order to complete the analysis of the cubic lattice, 
an investigation was also conducted into the behavior of the average latency 
period as a function of the parameters that governs the intensity and the time 
delay to mount a specific immune response. 

This paper is organized as follows: In section [2j the model is presented together 
with a discussion of the modifications introduced in order to deal with the 
triangular and cubic lattices. In section [31 the obtained results are presented 
and discussed, and the conclusions are presented in section HJ 



2 The Cellular Automata Model 



In the CA model [T] a target cell (T cells or monocytes) is associated to each 
site of a square lattice. Each cell is represented by a four- state automaton 
corresponding to different states of the cell: healthy, infected-A, infected-B 
and dead. The infected-A state corresponds to a virus-producing cell, which 
promotes the spread of the infection during a time interval r without any 
suppression, r corresponds to the period of time necessary for the immune 
system to mount a specific immune response to a given virus strain. A virus- 
producing cell that was already recognized by the immune system is then 
represented by the infected-B state. Although infected-B cells are less effec- 
tive than infected-A in contaminating healthy cells, a high concentration in 
a given region may disseminate the infection throughout the nearby unin- 
fected cells. The infected-B cells eliminated by the cytotoxic mechanisms of 
the specific immune response or the cytophatic effects of virus replication are 
represented by the dead state. The high ability of the immune system to re- 
plenish depleted cells is taken into account by converting a fraction p rep \ of 
the dead cells into new cells. Such new cells mimic the flux of incoming cells 
from other compartments. Most of the new incoming cells would be healthy, 
but a (small) fraction of them (p- m f ec ) would be infected-A, corresponding to 
virus-producing cells coming from other compartments. 

Simulations were carried on a L x L square lattice starting from an initial 
configuration composed of healthy cells, with small fraction, phiv> of infected- 
A cells randomly distributed among the healthy cells, representing the initial 
contamination by HIV. Periodic boundary conditions were considered and each 
time step corresponded to the parallel updating of the entire lattice according 
to the following rules: 

1) A healthy cell becomes an infected-A if it has at least one infected-A cell 
amongst its nearest neighbors, or at least R = 4 neighbors in infected-B 
state. Otherwise it remains healthy. 

2) An infected-A cell remains in this state during r time steps, after which it 
becomes infected-B. 

3) An infected-B cell becomes a dead cell in the next time step. 

4) A dead cell is replaced by a healthy cell with probability (1 — Pi n fec)Pre P i, 
or by an infected-A cell with probability p rep \ p- m f e c, otherwise it continues in 
the dead state with probability (1 — p rep i). 

Some of the parameter values chosen are based on experimental findings. For 
example, one in 10 4 to 10 5 cells in the peripheral blood of infected patients 
expresses viral proteins [16] during the latency period, and pinfec values are 
brought about by this finding. The probability of replacing dead cells by new 
incoming cells (p rep \) could vary between zero and one and may also vary 



from patient to patient. In the original model it was assumed that p rep \ = 
0.99 corresponds to a high replenishment capacity of the immune system, 
which in this case is not affected by the HIV infection, and is contrary to 
certain assumptions cited in the literature. Such a replenishment concerns to 
the population of cells that was recruited or committed do participate in the 
immune response. Finally, the small fraction phiv of infected cells randomly 
distributed in the initial configuration was chosen according to experimental 
findings, which indicates that one in 10 2 to 10 3 cells harbor viral DNA during 
primary infection [TT|[T8] . 



According to the updated rules adopted since the beginning, infection is dis- 
seminated in a deterministic manner, driven by cell-cell contact, but eventu- 
ally newly infected cells are introduced to the system due to stochastic rule 4. 
Rule 4 accounts for all mechanisms governing the (re) infection process, such 
as the presence of quiescent infected cells or infected cells coming to other 
compartments. 



In order to investigate the role of lattice symmetry and its dimensionality on 
CA model dynamics, studies were conducted for the models on the triangular 
(2D) and cubic (3D) lattices, to compare results with those previously obtained 
using the square lattice. For such cases, few changes should be made to the 
rules regarding the neighborhood, as described above: 



Triangular lattice: 

Each site possessed six neighbors, and R = 3 was assumed as the minimum 

number of infected-B neighboring cells required to disseminate the infection, 

i.e., following the same criteria adopted for the square lattice, this number 

would correspond to half the number of neighbor sites. 

Cubic lattice: 

Each site possessed twenty six neighbors, hence R = 13 was assumed as 

the minimum number of infected-B neighboring cells required to to infect a 

healthy cell. 



To test the robustness of the model's dynamic behavior regarding the varia- 
tion of stochastic parameters, the position and fraction of the infected cells at 
primary infection peak were measured, together with the duration of the la- 
tency period when varying respectively, phiv and pinfec, over wide intervals. The 
sensitivity of the latency period was also measured in the three-dimensional 
model, when changing the deterministic parameters R and r. These results 
are presented and discussed throughout the following section. 



3 Results and discussion 



In order to analyze the effects of the changes in the symmetry of the lattice, 
an average three-stage pattern was obtained that emerged from the CA model 
defined on the triangular and square lattices for L = 900 and the same set 
of parameters used in the original model pQ: phiv = 0.05, r = 4, p rep \ = 0.99, 
Pinfec — 10~ 5 . The results shown in figured] were an average of 1000 samples 
(initial configurations) corresponding to different individuals. Hereafter the 
error bars shown in the plots of all figures indicate one-standard deviation. 
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Fig. 1. Time evolution of the fractions of healthy (squares), infected- A + B (circle) 
and dead (triangles) cells. Solid symbols correspond to the results of the triangular 
lattice model (R = 3), while open symbols refer to the square lattice model (R = 4). 

As can be observed, there are no major quantitative changes in the patterns 
obtained for the two lattices (square and triangular), indicating that the three- 
stage dynamics is quite independent of the lattice symmetry. Figure [2] com- 
pares the results obtained from the model defined on the cubic and square lat- 
tices. Except for the changes in L and R, all parameters were those adopted 
in ref [TJ. It is worth observing that when the results for different dimen- 
sions were compared, significant timescale changes were observed. However, 
qualitative behavior remained the same in both cases. The differences of the 
three-dimensional results in relation to those of the two-dimensional were: a 
one week shift in the peak of the primary infection phase, and a reduction 
in the duration of the latency period. The duration of the latency period is 
estimated for each patient by calculating the time necessary for the fraction of 
healthy cells to achieve a threshold of ~ 30% after the primary infection. In 
this manner it was possible to estimate the point that corresponds to the onset 
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Fig. 2. Time evolution of the fractions of healthy (squares), infected A + B (circle) 
and dead (triangles) cells. Solid symbols correspond to the results obtained for the 
cubic lattice with L = 300 and R = 13 while open symbols correspond to the square 
lattice using L = 900 and R = 4. The parameters are those adopted to obtain Figure 
1. 

of AIDS [2]. The average length of the latency period (in over 1000 samples) 
was reduced from 7 (2-D) to 2 years(3-D). A decrease in the error bars mag- 
nitude was also observed. This result is easily understood when considering 
the number of neighbors in the different lattices. While the cubic lattice pos- 
sessed 26 nearest neighbors, the square and triangular lattices possessed 8 and 
6 neighbors, respectively. Hence, the large number of neighbors in the 3D case 
will favor a much faster dissemination of the infection, when compared to the 
2D lattices, thus reducing the latency period. When considering the lymphoid 
tissue as a fractal object, according to the results of this study, the fractal 
dimension should be closer to two, thus validating the approximation of a two 
dimensional lattice adopted in the original model [1] . The robustness of the 
dynamics of the CA model during primary infection was also investigated for 
variations of the initial concentration of infected cell (pmv) f° r 2D and 3D ver- 
sions of the model. Figure [3] and H] respectively, on a logarithmic scale, reveal 
the behavior of two quantities that characterize the primary infection phase: 
the average maximum concentration of infected cells < D{ > and its position 
< t > (in time steps) as a function of phiv- The small error bars indicate the 
robustness of the results irrespective of the lattice dimensionality. Power-law 
behavior was observed, indicating < Di >~ {pw\\i) a and < t >~ (phw/) 13 with 
exponents a = 0.48 ± 0.01 for two-dimensional square and triangular lattices 
and a = 0.29 ±0.01 for the cubic. These values of a are related to the average 
distance between first-nearest-neighbors infected-A cells in the initial config- 
uration, which by it turn scales with 1/D. Each infected-A cell of the initial 
configuration produces a pulse of infected cells, of width (r + 1) propagat- 
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Fig. 3. Average maximum fraction of infected cells during primary infection as 
a function of £>hiv- The different plots correspond to: cubic (0)> square (□) and 
triangular (A) lattices. The results were an average of over 1000 samples the error 
bars being of the order of the symbol size. Dashed lines indicate the linear fitting. 



ing in all directions. Whenever such average distance is less than or equal to 
(2r + 1), the independent pulses achieve a maximum coverage of the lattice, 
which corresponds to a maximum of < D L >. The values obtained for (3 are 
similar for the two-dimensional lattices: (3 = 0.44 ± 0.01, but a decrease was 
recorded of 0.30 ± 0.02 for the cubic. In other words, the primary infection 
peak becomes broader with the decrease of height by one order of magni- 
tude while phiv varies over two decades for two-dimensional lattices and three 
decades for three-dimensional. Other authors have claimed that there exists a 
lack of robustness of the CA model when describing primary infection |19j . 
However, the results of this study assert the robustness of the results obtained 
with respect to variations of the minimum concentration of infected cells nec- 
essary to launch the infection process. The claim was based on an unproven 
assumption of exponential growth of the viremia titer in the very beginning 
of primary infection, and on speculative values for the minimum amount of 
virus necessary to launch the contamination of an individual, inferred from a 
few clinical cases. This study illustrates that the same dynamics are launched 
for a wide range of the initial concentration of infected cells. 
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Fig. 4. Average position of the maximum fraction of infected cells achieved during 
primary infection (< t >) as a function of Phiv- Symbols represent the same lattices 
as described in Figure [3] and the calculation of averages and error bars also follows 
the same procedures as in the previous figure. Dashed lines indicate the linear fitting. 

The average latency period was also investigated as a function of the proba- 
bility of newly infected cells entering the system (p rep | x p; n fec) - Since the hy- 
pothesis was adopted that the function of the immune system is not affected 
by the infection, and the replenishment of healthy cells by the bone marrow is 
maintained constant (p reg = 0.99), we have plotted the average latency period 
solely as a function of pi n f ec , as shown in figure [51 The behavior of the average 
latency period is similar for both the triangular and square lattices, but is 
different in the case of the cubic lattice. The power law exponent 7 associated 
to the average latency period is 0.40 ± 0.01 for square and triangular lattices, 
and 0.25 ± 0.01 for the cubic. 

As previously mentioned, when the same set of parameters (except of L and 
R) are taken and the dynamical behavior exhibited by the two and three- 
dimensional models is compared, a significant reduction of the average latency 
period can be observed. Recently, Solovey et al. [20] investigated how the varia- 
tions on the deterministic parameters change the dynamics of the model in the 
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square lattice. To complement the analysis concerning the differences between 
the two and three dimensional models, this study investigated the behavior 
of the average latency period in the three-dimensional model when varying R 
and r. The results in figure [(^a) indicate that, as in the case of the square 
lattice [20] , the average latency period is not significantly altered due to vari- 
ations in the time-lag parameter r. However, when parameter R varies from 
9 to 18 an increase in the latency period is observed to the order of 40 weeks, 
as displayed in figure [6](b). Such an increase corresponds to a decrease in the 
probability of the occurrence of infected cells and concomitant compact struc- 
tures leading to long latency periods, as would be expected. It can be noticed 
that both plots of figure [6] exhibit large fluctuations. Such fluctuations are 
related to the intrinsic width of the latency period distribution and not with 
the number of simulations. In fact, fluctuations are governed by the spatial- 
and time- dependent probability of the appearance of compact structures, and 
the latency period is related to the average distance of such structures, which 
becomes proportional to the inverse of the cubic ratio of the concentration p 
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Fig. 5. Average clinical latency period as a function of pi n fec- The different plots 
correspond to: cubic (0)> square (□) and triangular (A) lattices. Results were an 
average of over 1000 samples. Dashed lines indicate the linear fitting. 
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Fig. 6. Average clinical latency period for the cubic lattice model: (a) as function 
of t; (b) as a function of R. Results were averaged over 50 samples, using L = 300. 

4 Conclusions 



This paper has analyzed the robustness of the cellular automata model pro- 
posed in reference [I] to describe the HIV infection with respect to changes 
in lattice symmetry and dimensionality and in variations of the stochastic 
parameters of the model. It has been shown that during the entire course of 
the infection, the time evolution of the fractions of healthy, infected and dead 
cells, are quantitatively reproduced irrespective of the lattice symmetry in two 
dimensions. However, the average latency period was found to be greatly re- 
duced when the lattice dimension changes from two to three, while the primary 
infection scale and the overall qualitative behavior remains almost unchanged. 
These results indicate that interactions within the lymph nodes occur on an 
effective "surface" with a fractal dimension close to two. This therefore ex- 
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plains the success obtained with the two-dimensional approximation used in 
the original model. 

The robustness of the three-stage pattern was also observed when varying the 
stochastic parameters associated to the minimum amount of infected cells in 
the initial configuration necessary to launch the infectious process (phiv) an d 
to the mechanisms responsible for either incoming newly infected cells or acti- 
vating quiescent infected cells (jOinfec)- The results indicate quite a remarkable 
robustness of the three-scale dynamic pattern following power-law behaviors 
over three orders of magnitude of Phiv- The observed width increase together 
with the decrease of the primary infection peak, indicates that the lower the 
initial concentration of infected cells, so the spread of the infection within the 
tissue is slower and weaker, as would be expected. 

The behavior of the average latency period as a function of pi n fec was also 
investigated for the three different lattices (square, triangular and cubic). The 
results also exhibit power-law behaviors over four (2D) to five (3D) decades: 
the lower the Pi n fec the greater the latency period, since the probability of form- 
ing structures associated to the rapid evolution of the onset of AIDS is reduced. 
In all cases studied, it was observed that the lasting latency period depends 
on the formation of spatial structures. Once these structures are formed they 
spreads all over the lattice, thus compromising tissue and trapping healthy 
cells, as observed in the square lattice model [T] . Such structures could be as- 
sociated to aggregates of infected cells observed in vitro experiments, namely 
syncytia, and in the lymph nodes of infected patients. Therefore, by means of 
the robustness analysis, the hypothesis can be validated that the formation of 
such structures may be responsible for the persistence of the virus in the sys- 
tem after primary infection [3"f2"T| . The appearance of such aggregates, which 
naturally emerge from the dynamics of the system in the automata cellular 
model, suggest that a deeper biological investigation should be performed in 
order to confirm the existence of these structures in vivo and their role in the 
HIV infection. 
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